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Self-consistent solutions of Hedin’s diagrammatic theory equations (HE) for the two-site Hub¬ 
bard Model (HM) have been studied. They have been found for three-point vertices of increasing 
complexity (T = 1 (GW approximation), Fi from the first order perturbation theory, and exact 
vertex Fb). The comparison is being made when an additional quasiparticle (QP) approximation 
for the Green function is applied during the self-consistent iterative solving of HE and when QP 
approximation is not applied. The results obtained with the exact vertex are directly related to the 
presently open question - which approximation is more advantageous for future implementations - 
GW-I-DMFT or QPGW+DMFT. It is shown that in the regime of strong correlations only orig¬ 
inally proposed GW-I-DMFT scheme is able to provide reliable results. Vertex corrections based 
on Perturbation Theory (PT) systematically improve the GW results when full self-consistency is 
applied. The application of the QP self-consistency combined with PT vertex corrections shows 
similar problems to the case when the exact vertex is applied combined with QP sc. The analysis 
of Ward Identity violation is performed for all studied in this work approximations and its relation 
to the general accuracy of the schemes used is provided. 
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I. INTRODUCTION 

One of the challenges for the computational theorists 
working in the solid state electronic structure field is 
the robust implementation of the so called GW-I-DMFT 
method (combination of GW approximation (G- Green’s 
function, W - screened interaction) and Dynamical 
Mean-Field Theory). The scheme was originally pro¬ 
posed by Sun and Kotliar— and in slightly different (but 
probably more commonly known) form by Biermann et 
alj^ The basic idea of the approach is to separate all ac¬ 
tive space of the basis set into ’’weakly correlated” part 
for which GW approximation is supposed to work well 
and ’’strongly correlated” part for which one sums up 
many diagrams (to infinite order if one uses the non- 
perturbative DMFT solver as for example the Quan¬ 
tum Monte Carlo - QMC). In order for the method to 
be well defined, everything should be done till full self- 
consistency (sc), including the iterations of the GW part 
itself and also the ” internal” iterations of DMFT part to 
ensure that the solution of the impurity problem repro¬ 
duces the same G and W in strongly correlated subspace 
as the ones in the same subspace projected from the GW 
part. 

Despite its obvious appeal GW-I-DMFT has made only 
slow progress during more than decade since its first ap¬ 
pearance in the above mentioned papers. The reason was 
not just because scGW is quite demanding computation¬ 
ally but mostly because one has to satisfy the impurity sc 
condition not only for G (as in the LDA-I-DMFT method 
- a combination of Local Density Approximation in Den¬ 
sity Functional Theory and DMFT) but also for the W, 
which seems to be not an easy task. To the best of my 
knowledge there were only ” one-shot” type calculations 
for real materials^ where GW iterations were neglected 
altogether and DMFT self-consistency was imposed only 


on G, whereas W was fixed at LDA level and correspond¬ 
ingly the U was considered as an external parameter 
(calculated in constrained random phase approximation 
- cRPA). Such an implementation, of course, is a definite 
step towards the full GW-I-DMFT scheme, but still one 
cannot say that there was the summation of all ’’corre¬ 
lated” graphs in it which would require that the W in GW 
part and Wimp in the impurity (DMFT) part were the 
same. Instead they were totally decoupled which makes 
it unclear of what kind of diagrams from DMFT part 
were actually added to GW part. Nevertheless, together 
with development of QMG solvers capable to handle dy¬ 
namical interactions^!^ the hope is growing that the full 
GW-I-DMFT scheme will eventually be implemented. 

There are some subtleties about the GW part as well. 
In its full sc implementation the method is very time con¬ 
suming which in part has prevented its applications for 
the real materials. But recently a very efficient imple¬ 
mentation of it was published^ where the most computa¬ 
tionally demanding parts (calculation of the polarizabil¬ 
ity P and the self energy E) are performed in real space 
and Matsubara’s time. As a result, it became possible 
to successfully apply scGW to the actinides Pu and Am 
and (earlier) to simple sp-materials^^ However for the ma¬ 
jority of solids, scGW produces worse spectra than fast 
” one-shot” GW and this is the other reason why scGW 
is not popular. The reason for the failure of scGW with 
spectra may be traced as an extremely non symmetri¬ 
cal ’’dressing” of the ’’initial” Green’s function with self 
energy insertions of GW-only form in the course of self- 
consistency iterations and neglecting by vertex correc¬ 
tions. 

The origin of the problem with the scGW method 
can also be formulated in terms of the absence of the 
Z-factor cancelation;^ which again happens because we 
neglect by the vertex corrections. To resolve this prob- 
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lem Kotani and Schilfgaarde^ devised a beautiful trick 
of doing yet another approximation. They used quasi¬ 
particle (QP) part of Green’s function Gqp (instead of 
full Green’s function G calculated from Dyson’s equation 
(DE)) to calculate P and E on every iteration till self- 
consistency. The trick is, that the errors from the above 
two approximations (using Gqp instead of full G and ne¬ 
glecting by vertex corrections) mostly cancel each other 
out and as a result the QPscGW (self-consistent quasi¬ 
particle GW) method usually gives much better spectra 
than full scGW. The important fact is that QPscGW not 
just slightly improves the one-shot GW description of sp- 
materials (which are good enough already in the one-shot 
GW), but often gives reasonable results for the materials 
with d- or f-electrons too,— and the method doesn’t rely 
on a particular starting point. It is totally self-consistent. 

Quite naturally, the success of QPscGW with spec¬ 
tral properties (as compared to the full scGW) has 
ignited the ideas to formulate another approach - 
QPscGW-l-DMFT— iii, where one supposedly adds 
DMFT corrections to P and E as in GW-I-DMFT but 
uses QPscGW for the "big” iterations in the ’’weakly 
correlated” part. 

As it appears we have now two schemes proposed: 
GW-bDMFT and QPGW-bDMFT. In this work I am 
doing an attempt to ’’estimate” what to expect from 
the future implementation of the schemes. The analy¬ 
sis strongly depends on the ” separability” of weakly and 
strongly correlated parts. I will assume here that they are 
perfectly separable and the correlations in the GW part 
are really weak. In this case, Z-factor in the GW part 
is close to 1, so that this part is equally well described 
in both the GW and the QPGW approximations. The 
difference correspondingly comes only from the DMFT 
part. For this part I assume that we are able to solve 
it exactly which means in particular that both sc condi¬ 
tions (for G and W) are satisfied. For GW-I-DMFT, it 
means G = Gimp, W = Wimp, and for QPGW-bDMFT, 
it means Gqp = Gimp, Wqp = Wimp, where following 
the arguments of work Q we have Gqp ~ The exact 
solution of the impurity problem also means that corre¬ 
sponding self energies can be written in their exact dia¬ 
grammatic forms: E = —GTW and Eqp = —GqpTWqp 
with r being the three-point vertex function which is in¬ 
cluded exactly. But following again the arguments in 
Ref.Q, it becomes clear that in QPGW-I-DMFT case, i.e. 
when we include exact vertex and continue to apply QP 
approximation for G, we will have a problem, because the 
factor 1/Z appears from the vertex and it doesn’t can¬ 
cel with the Z-factor in G as it happens in GW-I-DMFT. 
Basically it means the violation of Ward Identity (WI) 
in the QPGW-I-DMFT scheme. Thus, from this point 
of view, QPGW-I-DMFT is highly problematic and its 
difficulties should grow up when the correlation strength 
grows up, because the l/Z-factor increases. 

The above simple argument against using the 
QPGW-I-DMFT scheme for strongly correlated materi¬ 
als needs some numerical support and I will provide it 


here using the exactly solvable two-site Hubbard model. 
For this model I calculate exactly the three-point ver¬ 
tex function and use it to calculate self-consistently the 
Green functions G and Gqp following two slightly differ¬ 
ent sc schemes (the scheme on the left side is basically 
the Hedin’s sc equationsiSi but with the 3-point vertex 
precalculated exactly) 

P — GY pxactG Pqp — GqpY pxactGqp 

W = u + UPW Wqp = 17 + UPqpWqp 
E = —GTExactW Z,QP = —GqpTExactWQp ^ ^ 
G = Go + GoEG Gqp = ^j^S^, 

where U is the bare interaction in the Hubbard model 
and Go if Green’s function in Hartree approximation. In 
o for QP case I formally represented the quasiparticle 
approximation for G by simply dividing the full Green 
function from Dyson’s equation by the Z-factor. This 
is for brevity. In fact I use the algorithm described in 
Ref.@ to construct Gqp. 

It is obvious that at self-consistency the scheme on the 
left side (I will call it GP eW) is equivalent to the ex¬ 
act solution of GW+DMFT equations whereas the right 
hand scheme (I will call it QPGT eW) is equivalent to 
the exact solution of the QPGW+DMFT equations. It 
is also obvious that the GT eW scheme is exact by con¬ 
struction. I use it in this work to check the numerical 
accuracy of 3-point vertex evaluation. The QPGV eW 
scheme is approximate and below I will explore its accu¬ 
racy in different regimes of correlation strength for the 2- 
site Hubbard model. Also I will directly relate the prob¬ 
lems of the QPGT eW scheme with the degree of the WI 
violation. 

Another goal of the present work was to explore the 
possibility of combining the GW and the QPGW meth¬ 
ods with perturbative calculation of the 3-point vertex 
function. To this end I will use again the schemes similar 
to Eq.([T]) but with P expanded to the first order in W (Pi) 
instead of exact P^;. I also will show how the two corre¬ 
sponding perturbation theory based schemes (GPiIT and 
QPGTiW) behave in the different regimes of parameters 
of the Hubbard model. 

This paper begins with the formal presentation of the 
two-site Hubbard model and the formulae used in the 
calculations (Section |H| . In Section |HI] the results are 
presented and discussed. Finally in Section |TV] the con¬ 
clusions are given and the future plans are outlined. 

II. TWO-SITE HUBBARD MODEL 

The Hamiltonian of two-site Hubbard model as it is 
used in this work is the following: 

H = —t + — '^^'^^cf^cf^iCia-’Cicr, ( 2 ) 

i <7<7' 
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where t and U are the standard parameters of the Hub¬ 
bard model, c and are the destruction and creation 
operators correspondingly, indexes i and j belong to the 
sites (1 or 2), and a, a' are the spin indexes. 

In this section, all equations which are used in the anal¬ 
ysis of the 2-site Hubbard model are collected for refer¬ 
ences. First I provide the energies and eigen vectors of 
the exact many-body states of the model for different oc¬ 
cupancies. Then, the exact expressions for Green func¬ 
tion and density-density correlation function are given. 
From them the exact self energy, polarizability, dielec¬ 
tric function, and screened interaction can be calculated 
using standard formulae. Next, two subsections provide 
the formulae for the exact 3-point vertex function and 
how it is used to evaluate the corrections to the polar¬ 
izability and to the self energy in sc scheme O- Per¬ 
turbation theory based equations for the 3-point vertex 
function are given next. Then I provide full and reduced 
(long-wave and long-wave-|-static limits) expressions for 
the Ward Identity which are used later in this paper. 
Finally the formulae relevant to the evaluation of cur¬ 
rent three-point vertex function are given. For brevity, 
the derivations of the formulae are not provided, or only 
sketch of derivation is given. In this work, the finite tem¬ 
perature framework is used so that in all equations be¬ 
low the time-frequency arguments are Matsubara’s time 
(r) and Matsubara’s frequencies (w for fermion frequency 
and v for boson frequency). 


A. Many-body states 


for N=2 

£^2 _ C/-c.,j-,2 _ |H1001)-|0110) |1100) + |0011) 

1 2 > 1 y/a{c-U) 

El = 0;^-! = 11010) 

El = 0 ;^'i = joiOl) 

F;|=0;^1 = ^(11001)+ 10110)) (5) 

iff = C/;'ki = ^(11100)-10011)) 

t^2 _ U+c. ,t, 2 _ ^^|U)0]+d0n+ |1100)+|0011) 

® 2 ’ 6 ^6(c-h(7) ^ ’ 


with a = ^2 + b = ^2+ and c = 

y/m?TTP. 

For N=3 

if? =t/-t;vI/3= ^(11011)-11110)) 
iff = t/-<;vI/3 = ^(10111)-11101)) 
iff = t/ + <;vl/i = ^(11011)+ 11110)) 

El = U + t- = ^(10111) + IllOl)), 

and for N=4 


E* = 2U-'i’f = \nn). (7) 


B. The Partition function, Green’s function and 
self energy 


For convenience, first I introduce ’’shifted” many-body 
energies E!^ = E^ — ^N, with /r being the chemical 
potential. Then I renormalize them, defying the minimal 
among them and subtracting it if"^ = E'J^ — E'^^^. 
This also factorizes the partition function: 


In order to represent the many-body states of the 
model it is convenient to introduce basis vectors \abcd) 
where all entries are equal to 0 or 1 in accordance with 
the occupancies of corresponding one-electron states. In 
this work first two one-electron occupancies (a and b) 
correspond to spin-up and spin-down one-electron states 
of the first site, and the third and fourth (c and d) corre¬ 
spond to the second site. Many-body energies and states 
below have two indexes: the upper one corresponds to 
the full occupancy (0,1,2,3, or 4) of the system, and the 
lower one distinguishes the states within the same occu¬ 
pancy. 

For the full occupancy N equal zero, we have corre¬ 
spondingly: 


F;? = = 10000), (3) 


where /3 is the inverse temperature. 

It is clear that now in every Gibbs average one can use 
E"\ Z' instead of E'\Z which is numerically more stable 
(big numbers have been subtracted). 

For exact Green’s function, one obtains through the 
standard spectral decomposition: 


= E ^ 

N n^N m^N+1 ^ mN+1 ^ nN 


X (4/):|c..K+^)(vi- 


N+l\ 




N\ 


(9) 


for N=1 


The exact self energy is obtained by inversion of the 
Dyson’s equation: 


F;i = -t;4'i = ^(|1000)+ 10010)) 
F;i = -t;4'^ = ^(10100)+ 10001)) 
^3 =i; ^3 = ^(11000)-10010)) 
^4 = i;'1^4 = ^^(10100)-10001)), 


S);(a;) = Go;j(cc)-G-.'-(ca), (10) 

(4) where Green’s function in Hartree approximation 
= (iw + /r - Upi)5ij + t(l - 5ij) is used {pi 
is the occupancy (’’density”) of the site i). 
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C. Response function, polarizability, dielectric and the screened interaction 

function, and W 

The exact two-point density-density correlation func- = USij + U^Rff{v). (15) 

tion is also obtained through the spectral decomposition 


dd 




.-PK 


— e 


-/3B" 


Z' ^ ^ ^ 11/ - 

N neN meN 




T^fi _ 77’// 

^nN ^mN 

-rN\ 


( 11 ) 


with density operator Pi = H 1® convenient to 

define also the density-density response function 


The response function, the dielectric function, the po¬ 
larizability, and the screened interaction calculated using 
the formulae (IT^ - (I15|) from the exact correlation function 
are by construction exact and can be compared with the 
corresponding quantities obtained using the PT. 


Rffiu) = (12) 

After that one can find the density-density dielectric 
function 

e.y^"(^)=^..+tl<(^), (13) 

the density-density polarizability 

^"(^)=E^''(^)4'(^)’ (14) 


D. Exact 3-point vertex function in density channel 

To find the exact 3-point density vertex function I first 
calculate the following three-point correlation function: 

= {cta{T)c+^{T')pk{0)) with t,t' being Mat- 
subara’s times and ’d’ meaning ’density’. In the site- 
frequency domain its spectral decomposition reads as the 
following: 


J 


Xy’fc(w;^) = ^E{ 

N 

-E E 






pSAf raeV+l 

+ E E 


mSN+l 

p€N+l ra^N ~*(^ — v) + Ep — E!^ 


IV + -tjn ~ tL-m 

^-BE‘ 


lU 


+ E''N - E''N 




= -/3B" 


-f 




+ E''N - 


ILU 


+ E E 

p£N-\-l m£N 








n£N 


iuj + E'l^ - E'4+^ 


}• (16) 


After that the three-point density response function is 
calculated 



5(j)k{E) 


pS^oG'^j{uj)pk + xlfki^'^ ^)’ 

(17) 


4>k{E). The ’’screened” 3-point density vertex function 
Xijk^'^-’ defined as the functional derivative of the 
inverse Green’s function with respect to the external per¬ 
turbing field (jjk (y) and correspondingly can be related to 
the above defined 3-point response function 


where I have indicated that the three-point response 
function is defined as the functional derivative of Green’s 
function with respect to the external perturbing field 
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_6G^l^ 

Hk{v) 

Y, z.)G;;-(a; - i.). (18) 

It 


Finally the three-point vertex function entering the 
Eq.Q (’’bare” three-point vertex) is defined as the func¬ 
tional derivative of the inverse Green’s function with re¬ 
spect to the total field <hfc(:^) (perturbing external plus 
induced internal) and is related to the ’’screened” ver¬ 
tex through the density-density dielectric matrix (this 
equation is the density-density part of the more general 
equation (IdTll l 


<5G:, '^(u;) 


ncr,d/ , \ _ 

- -- 

_ \ cr,cZ 


I 




(19) 


E. Vertex corrected polarizability and self energy 

The vertex-corrected density-density polarizability and 
the self energy entering the equations o are calculated 
as the following 

(^) = ^ E E E i.)GUu; - i^), 

^ UJ <7 kl 

( 20 ) 

and 

E E y)Wu{y). 

^ V kl 

( 21 ) 

The formulae (I20|) and (EU are used also when the 
three-point vertex is obtained within the perturbation 
theory. 


F. 3-point vertex function from Perturbation 
Theory in density channel 

The first order (in W) term of the perturbation theory 
for the three-point density vertex function is 

KfkiGT') = -lF,,(r' - r)G(;(r)G^,■(-/). (22) 

In the frequency representation it can be conveniently 
evaluated as the following 

(..; ^) = - / rfrlT,,(r)i ^ e—Grfe(a;)G^,(a; - u). 

UJ 

( 23 ) 


G. Ward Identities 

In order to write down the Ward Identities one needs to 
specify the current operator. It is convenient to introduce 
it through the substitution for the kinetic part of the 
Hamiltonian 


H' = -tY Y E E 'Grr'Ga, 

i(T i aa' 

(24) 

where Ai is the vector potential and the following con¬ 
vention for the sites was adopted: 1 = 2;2 = 1. With 
the above definition the current operator for the model 
reads as the following 

* CT 

with Ji = —Ji- From the equation of motion for the 
density operator one calculates 

^7 = (26) 

which is the continuity equation for the model. Relat¬ 
ing the ’’screened” 3-point current vertex function with 
the corresponding 3-point response function (similarly to 
the density case and in the space-time coordinates for 
brevity; ’c’ goes for ’current’)) 


G'"(14)7'"’7453)G‘"(52) = G'"(12)J(3) 

+ (c71)4(2)J(3)), (27) 

one explicitly calculates the time derivatives and using 
also the continuity equation (1261) one obtains 

G'^(14){E^M(453) + i7'-'7453)}G'"(52) 

= E(c71)c+(2)p( 3)) +i(c7l)c+(2)J(3)) 

= G'"(12){d(13)-<5(23)}, (28) 

which is equivalent to 


E-j,<^.<i(123) -p f7‘^’‘’(I23) = G '‘^(I2)|(5(23) - d(I3)|. 

(29) 


In the site-frequency representation (|2^ reads as the 
following 


(w;i^) + 75fc(w;j^)} 


“ Gij '^i^)Sjk — Gj^j — v)Sik. ( 30 ) 
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The equation (1301) can be simplified by removing the 
Hartree-Fock contribution on the both sides of it. For the 
vertices, one obtains in the Hartree-Fock approximation 

(31) 

and 

- SjkY (32) 

For the schemes with full sc (without the QP approx¬ 
imation) the removing of the Hartree-Fock contribution 
also on the right side of (1501) through Dyson’s equation 
gives 


v) = E)'/(w - v) - S7'^(w), (37) 

k 

and 

iv ^ (w; v) = 0. (38) 

k 

From (1551) one can see that in order to satisfy the long¬ 
wave limit of the WI in the QP-based approximations 
one has to neglect by vertex corrections altogether. 

The long-waved-static limit of the WI consists in taking 
the limit (;/—>• 0) in the equations (157)) and (1551) : 


= (33) 

where Ay means the "screened” vertex part beyond the 
Hartree-Fock approximation, and the E'’ is the correla¬ 
tion (frequency dependent) part of the self energy. For 
the QP-based schemes the Dyson equation is not satisfied 
and instead of (1551) one has 

*{'^Ayfy^(a;; v) + Ay'^.’(;(a;; z/)| 

= [h^/ (34) 

where the static effective Hamiltonians and 

•'J ‘'J 

(correspondingly in the Hartree-Fock and in the QP ap¬ 
proximation) were introduced: 


= = 
k 

and 


lim 




(39) 


(w; = 0) = 0. (40) 

k 

The arguments supporting the quasiparticle approxi¬ 
mation in Ref. @ are based on the long wave limit of the 
Ward Identity. It is clear from the above consideration 
that in the QPGW approximation (without the vertex 
corrections) the corresponding limit of the WI is satis¬ 
fied exactly. 


H. 3-point current vertex function 


= -t(I - Jy ) + + Sn, (35) 

with RA and Ef being the Hartree potential and the 
exchange part of the self energy correspondingly, and 

H?/ = ^ 

kl 

(36) 

where Zij is the renormalization factor, and the E^;(0) 
stands for the correlation part of the self energy at zero 
frequency. When the effect of correlations is neglected 
the renormalization factor becomes equal to I, and the 
correlation part of the self energy becomes zero, which 
means that in this case . 

The equations (1551) and (1551) are used in this work to 
evaluate the deviations from the full Ward Identity for 
different approximate methods. I also use two reduced 
forms of the WI in the present study: the long-wave limit 
and the long-wave-l-static limit of the WI. The long-wave 
limit of the WI for the two-site Hubbard model consists in 
the summation over the index k in the equations (1331) and 
(1341) which correspondingly become (the current vertex 
disappears after the summation) 


In order to apply the full WI one needs the ’’ screened” 
current vertex function y'^’(?(w; v). In this work both the 
exact and the PT-based vertex functions are used. The 
exact one is evaluated following the formulae dm-dni) 
with the replacement (p —J) in the Eg. dTBl) . The PT- 
based vertices are calculated following the scheme out¬ 
lined below. 

The "bare” vertex functions are related to the 
” screened” ones through the full (density-current) dielec¬ 
tric function via the following equation which is the gen¬ 
eralization of Eq. (ITOl) 

Tzy’fc (w; (‘^^ (41) 

J i 

where both I and J now run over the indices d (density) 
and c (current). The non-perturbed Hamiltonian has no 
vector potentials so that the full dielectric matrix has the 
form (considering the first index as the density, and the 
second as the current) 


e = 



^dc 


1 


(42) 
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and correspondingly its inverse 

^dd _ ^~^dd^dc \ 

0 1 j- 



(43) 


Thus, the bare current vertex can be evaluated as the 
following 




The subtraction of the Hartree-Fock contribution (l32l) 
leaves us with the expression 



A75=(..;z.)=Ar-:f(a;;^) 

- Y^lSiiSji + Ar^jf (w; 

Im 

(45) 

In the above expression the missing components are the 
non-trivial part of the current vertex Ar(^’°(w; v) and the 
density-current dielectric matrix e^fc(!^)- The first one 
is evaluated similar to the equation (1231) for the density 
vertex: 


FIG. 1: (Color online) Imaginary part of the on-site self en¬ 
ergy as a function of frequency for different average occupan¬ 
cies for U=l. 


(49) 

nN 

To evaluate the internal energy in the perturbation the¬ 
ory based methods the following formula is used 


^ U) 

X {Gf,(cc)G?^.(cc -v)- Gl{uj)Gl^{uj - i^)}. (46) 

The second one is defined by the density-current po¬ 
larizability P‘^'^ 

eg=(z.) = -UPffiu), (47) 

which in its turn is evaluated using the current vertex 

(HU): 

= zt f drP’'- Y. {Gy mg? (-r) - G-.(t)G-(- r)} 

(T 

+ ^ E E E ^)GU^ - U). (48) 

^ OJ cr kl 

I. Internal energy 

The exact internal energy is evaluated directly as the 
average value of the Hamiltonian. In the spectral repre¬ 
sentation it reads 


A E GJ.(t = /3) + ^ ^ ^ P^aP^a' 

i aa' 

-^EEsMGr.(r = ^) 

i a 

(50) 

^ ia LO 

which is based on the Galitskii-Migdal expression for 
the exchange-correlation energy. 

III. RESULTS 

The two-site Hubbard model is studied here with the 
value of parameter t fixed and equal to 1. The tempera¬ 
ture was also fixed at T=0.05t. Thus, only the parameter 
U was changing. The half-filling {{N) = 2) case has been 
considered. Such a choice for the occupancy has been 
guided mostly by the fact that when one steps aside from 
the half-filling the correlations in the model quickly be¬ 
come unmanageable for the QP-based approaches. It is 
seen in the Fig IT] where the imaginary part of the on¬ 
site self energy is plotted as a function of Matsubara’s 
frequency for the average occupancies 1, 2, and 3. 










FIG. 2: (Color online) Errors in Green’s function for U=0.5. 


FIG. 5: (Color online) Errors in Green’s function for U=3. 




Frequency 

FIG. 3: (Color online) Errors in Green’s function for U=l. 

FIG. 6: (Color online) Errors in Green’s function for U=4. 
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Frequency 

FIG. 4: (Color online) Errors in Green’s function for U=2. 



Frequency 


FIG. 7: (Color online) Errors in Green’s function for U=5. 
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The strong downturn in the function for the occupan¬ 
cies 1 and 3 at the small frequencies makes the lineariza¬ 
tion in the self energy pertinent to the QP approxima¬ 
tions highly inaccurate. At the half filling however, the 
model shows slow increasing in the correlation strength 
with the increasing of U and, correspondingly, is conve¬ 
nient for the comparative studies. The calculations have 
been performed for different values of U which were in¬ 
creased till the methods based on the QP approximation 
began to fail seriously. 

I compare the exact results with the results obtained 
with the GW, the GPiW, the GP^W, the QPGW, 
QPGPiW, and the QPGPbW methods, where Pi and 
Pe stand for the first order (in W) 3-point vertex and 
the exact 3-point vertex correspondingly. At the full self- 
consistency GPeW reproduces the exact results, so this 
approach was used basically only as a mean to adjust 
(by comparing with the exact result) the calculational 
parameters (such as the number of Matsubara’s frequen¬ 
cies included in the internal summations and the density 
of mesh on the interval [0 : /3] for the r-integrations. 
In addition to the Green function which serves in this 
work as basic representation quantity for the compar¬ 
isons, the analysis of the Ward Identities fulfillment has 
been performed and was used as an indicator of the accu¬ 
racy of the approximations. Besides, the internal energy 
has been evaluated and its accuracy was related to the 
errors in the calculated G and to the deviations from the 
Ward Identities. 

In the figures [2]1Z1 the absolute error in the calculated 
Green’s function is shown as a function of the Matsub¬ 
ara’s frequency for the different values of U. The com¬ 
parison is being made between the exact result and the 
results obtained with the approximate methods. The in¬ 
ternal energy as a function of U is presented in the FigH) 
The figure IHl shows the relative violation of the full WI 
(Eqsl33] and [SH) in different approximate methods as a 
function of U. The following form of the average devia¬ 
tion from the identity has been used 



u 

FIG. 8: (Color online) Internal energy as a function of U 
parameter. 


One can do a few observations from the results pre¬ 
sented. 

Most important for the present work is an observation 
following from the errors in Green’s function, that the 
QPGPeW method doesn’t show any noticeable improve¬ 
ment as compared to the other approximate methods 
studied. Opposite to that - with U increasing it quickly 
becomes the worst of the methods considered (especially 
for the higher frequencies). The observation perfectly 
correlates with the results obtained for the internal en¬ 
ergy and for the Ward Identities (especially with the re¬ 
sults for the full WI). It means that putting the exact 
vertex together with the QP self-consistency is not com¬ 
patible. 

Grouping the methods based on the QP-approximation 
in one group and the rest of the methods in another 
(they satisfy the DE), one can do another conclusion. 


T,tjk iLHS)^jkiuj; v) - {RHS)ijkiuj; v) 


V 1 


, (51) 


where {LHS) and {RHS) are the left- and the right-hand 
sides in (I33l34p correspondingly. The 3-point vertex per¬ 
tinent to the specific approximation (i.e. for instance the 
vertex Ay is zero in the GW and the QPGW) was used 
to evaluate the deviation from the WI in all approximate 
schemes, and the exact vertex was used in GPeW and 
QPGPeW. The summations over the frequencies w and 
V in 6ID were performed for |a;| <20 and 0 < e < 20. 

The figure (TU] shows the similar plot obtained with the 
long-wave limit of the WI (I37I38I) . The formula analo¬ 
gous to (ISTI) has been used but without the fc-summation. 
Finally in the figure [TT] the deviation from the long- 
wave-fstatic limit of the WI (I39l40p is shown. In this 
case, the summation over v has not been included. 
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FIG. 9: (Color online) Full Ward Identity average violation 
as a function of U parameter. 
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Namely for the frequencies larger than approximately 1.5 
the QP group has considerably larger errors in G than 
the group with the full self-consistency. This again cor¬ 
relates strongly with the internal energy graph and with 
the full WI. 

Further, in the methods satisfying the Dyson equation 
the increasing in the accuracy of the vertex (1 —>■ Fi —>■ 
F e) consistently improves the accuracy of G, finally mak¬ 
ing it exact when the exact vertex is being applied. Nat¬ 
urally, the same tendency can be noticed looking at the 
internal energy graph and at all WI graphs. 

Opposite to that, in the QP-based methods the situ¬ 
ation is more complicated. At lowest U (0.5 and 1) the 
QPGW is the best and the QPGF_eW is the worst at low 
frequencies with the reversed tendency at high frequen¬ 
cies. For U equal 2, 3, and 4 one can notice that now 
the QPGF£;W is the best and the QPGW is the worst at 
low frequencies with the reversed tendency at high fre¬ 
quencies. For U equal 5 the accuracy of the QPGF£;W 
approach deteriorates also at low frequencies and its fail¬ 
ing at higher frequencies becomes severe. The important 
point here is, that the above high-frequency tendencies 
in the QP-based methods correlate well with the tenden¬ 
cies in the internal energy as one can easily see. Namely, 
for U<1 the the QPGP^jW gives the best (among the 
QP-based schemes) internal energy and for U>1 it is the 
worst among them. This also correlates very well with 
the deviations from the full WI. 

For U=0.5 the QPGW approach produces the best G 
at lowest frequencies among all approximate methods. 
One can speculate, that this fact is similar to the well 
known fact, that in weakly correlated real materials the 
QPGW is probably the best approach for studying the 
spectra (which are defined by the low frequency behavior 
of Green’s function). This fact can be related with the 
exact fulfillment of the long-wave limit of the WI in the 
QPGW approach. 

An interesting observation comes from the figure for 


M 

S 

e 


O 

U 

H-l 


C 

O 

•H 

fd 

•H 

> 

( 1 ) 

Q 


0.7 

Exact, 

gQw, QPGW - 

' 

0.6 


GW . 

GF^^W . 

, 



QPGF^W 


0.5 


QPGrj,W . 


0.4 




0.3 




0.2 



. 

0.1 


. . 


0 


. . 




0 1 2 3 4 5 


U 


FIG. 10: (Color online) q ^ 0 Ward Identity average viola¬ 
tion as a function of U parameter. 
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FIG. 11: (Color online) q —>■ 0, ^ 0 Ward Identity average 

violation as a function of U parameter. 


U=I. As one can see there is a strict separation in the 
accuracy of G obtained within the QP-group and within 
the non QP-group. The reasons for this feature have not 
been studied in this work but it could be related to the 
above mentioned swapping between different tendencies 
in the QP-based methods at this U. 

The GW-based (without the QP approximation) meth¬ 
ods give numerically exact internal energy till 17 ~ 2 
which is not very surprising because they are conserving 
(in Baym-Kadanoff sensei^) and in the weakly correlated 
regime should produce accurate total energies. Generally 
it is not a good idea to use the QP-based methods to eval¬ 
uate the total energy because even if they originate from 
the same ik-functionali^ as the GW-based methods the 
Dyson equation (DE) is not satisfied in them anymore. 
As it was indicated above, the deviations from the exact 
Green’s function in the QP-based approaches are espe¬ 
cially noticeable at high frequencies which are important 
for the total energies evaluation. The results shown in 
FiglU clearly support this point of view. 

Internal energy obtained in QP-based approximations 
deviates significantly from the exact one even at small 
U. It is interesting that in the QPGW scheme the er¬ 
ror remains almost unchanged till U = 5 which prob¬ 
ably is accidental because at C/ > 3 the PT be¬ 
comes unreliable. The important difference between 
the GW-based and the QPGW-based methods is how¬ 
ever that the GFiW improves the internal energy as 
compared to the GW, whereas applying the vertex of 
improved accuracy (1 —>■ Fi —?> F^;) in the sequence 
QPGW—^-QPGFiW—^-QPGF£;W makes the results worse 
and worse which also tells us that one should not apply 
the vertex correction combined with the QP approxima¬ 
tion. 

As it is clear, the degree of violation of the WI corre¬ 
lates well with the errors in the Green functions at av¬ 
erage and high frequencies, which are responsible for the 
accuracy of the calculated total energies. And indeed. 
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the comparison of the Figures [8] and |9] tells us, that the 
accuracy of the total energy and the accuracy of the full 
WI fulfillment are closely related. Thus, the full WI can 
be useful as a measure of the accuracy of the calculated 
total energies. 

There was a hope, that the deviation from the long- 
wave+static limit of the WI correlates well with the er¬ 
rors in G at low frequencies. But for the two-site Hub¬ 
bard model this seems to be not the case (besides the 
above mentioned success of the QPGW approach at the 
smallest values of U). However, in real materials where 
the spectra obtained with the QPGW are generally no¬ 
ticeably better than the spectra obtained with the scGW 
the situation might be different. 

IV. CONCLUSIONS 

In this work the Hedin’s equations!^ for the two-site 
Hubbard model have been solved self-consistently with 
and without applying the quasiparticle approximation for 
the Green function. The study has been performed both 
when the exact three-point vertex function was used as an 
input and when the perturbative theory (in its zero and 
first orders in W) was used to evaluate the correspond¬ 
ing vertices self-consistently. The results of this work ob¬ 
tained with the exact vertex have direct impact on what 
one can expect from the future implementation of the sc 
GW-kDMFT and the sc QPGW-kDMFT schemes. As 
it has been shown here, only the GW+DMFT approach 
can be considered as useful approximation. However, as 
it was said in the Introduction, this work deals with an 
ideal situation, when the GW part and the DMFT part 
are perfectly separable and the subspace where GW is 
used is very weakly correlated (so that GW and QPGW 
give identical results for the weakly correlated subspace). 
In practice, it is not always the case. The correlations in 


the ” GW”-subspace might be noticeable. In such situa¬ 
tion, the QPGW might be superior (with respect to the 
GW) for the subspace not included in the DMFT part. 
The conclusions about the DMFT part obviously remain 
as before - one should not impose the QP approximation 
on G in the DMFT part. As it seems, in such circum¬ 
stances the preference should be given to the approach 
(GW-based or QPGW-based) depending on which sub¬ 
space (the ’’weakly” or the ’’strongly” correlated) is more 
important for the problem under consideration. However, 
on the fundamental level, such a situation should be re¬ 
solved either by the increasing of the subspace covered 
by the DMFT part or (which seems to be easier practi¬ 
cally) by including more diagrams beyond GW for the 
’’weakly” correlated subspace. 

It has been shown, that the methods with the PT- 
based vertices (when they are applicable) reveal similar 
tendencies (for example if one chooses between the QP 
self-consistency and the full sc) as the methods based 
on the exact vertices. Namely, when the correlation 
strength increases, both the PT-based and the exact 
vertices-based schemes begin to fail if one uses the QP 
self-consistency. Also of practical importance is the find¬ 
ing that the violation of the WI in any particular method 
correlates well with the general applicability of the given 
method. This can be an useful information if one sees 
to apply the schemes with vertex corrections to the real 
materials where the exact solutions are not available to 
serve as a judgement. 

The results of this work are of the methodological im¬ 
portance. Of course the two-site HM doesn’t cover all 
possible regimes of correlations which may happen in re¬ 
alistic materials. In order to cover a little bit more of the 
possible regimes of correlations the similar work on the 
homogeneous electron gas is now being performed (with 
the three-point vertex functions calculated within the PT 
only). 
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